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ABSTRACT 

We investigate the formation and subsequent co-evolution of galactic disks immersed in assembling 
live dark matter (DM) halos. Disk and halo structural components have been evolved from the 
cosmological initial conditions and represent the collapse of an isolated density perturbation. The 
baryonic component includes gas (which participates in star formation [SF]) and stars. The feedback 
from the stellar energy release onto the ISM has been implemented, and so is the mass, momentum 
and energy balance between the stellar and gaseous components. We find that (1) The growing 
triaxial halo figure tumbling is insignificant and the angular momentum (J) is channeled into the 
internal circulation, while the baryonic collapse is stopped by the centrifugal barrier; (2) Density 
response of the (disk) baryonic component is out of phase with the DM, thus diluting the inner halo 
flatness and washing out its prolateness over the time period of the disk growth; (3) The total J is 
neathly conserved, even in models accounting for feedback from stellar evolution; i.e., the baryons 
lose ~ 25% — 30% and the DM acquires ~ 2.5% — 3% of their original J; baryons and DM within 
the disk radius have the same J after the first ^ 3 Gyr; (4) The specific J for the DM is nearly 
constant, while that for baryons is decreasing with time; (5) Early stage of disk formation resembles 
the cat's cradle — a small amorphous disk fueled via radial string patterns; followed by growing oval 
disk whose shape varies with its orientation to the halo major axis, and so is the strength of a pair 
of grand-design arms developed in the disk; (6) The initially puffed up gas component in the disk 
thins when the SF rate drops below ^ 5 M© yr~^; (7) About 40%-60% of the baryons remain outside 
the SF region of the disk or in the halo in the form of a hot gas by the end of the simulations; (8) 
Obtained rotation curves appear to be flat and account for the observed disk/halo contributions (the 
halo virial mass is ~ 7 x 10^^ Mq at z = 0); (9) A range of bulge-dominated to bulgeless disks was 
obtained, depending on the energy thermalization parameter, esF, which characterizes the feedback 
from the stellar evolution — smaller esp leads to a larger and earlier bulge; Lower density threshold 
for SF leads to a smaller, thicker disk; gravitational softening in the gas has a substantial effect on 
various aspects of galaxy evolution and mimics a number of intrinsic processes within the interstellar 
medium; (10) The models are characterized by an extensive bar-forming activity within the central 
few kpc whose properties vary with the bar-halo orientation; (11) Nuclear bars, dynamically coupled 
and decoupled form in response to the gas inflow along the primary bars, as shown by Heller, Shlosman 
& Athanassoula. 

Subject headings: galaxies: evolution ~ galaxies: formation - galaxies: halos - galaxies; kinematics 
and dynamics - galaxies: structure - cosmology: dark matter 
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1. INTRODUCTION 

Hierarchical Cold Dark Matter (CDM) cosmology has 
been largely successful in establishing the paradigm of 
structure formation and evolution on spatial scales of 
~ 100 kpc, corresponding to individual halos, and be- 
yond (e.g., Cole et al. 2000; Steinmetz & Navarro 2002; 
Spergel et al. 2003). However, on smaller scales, a 
number of fundamental issues related to the origin of 
luminous galactic components and their structural pa- 
rameters remain unresolved because of our limited un- 
derstanding of dynamical processes, such as mass and 
angular momentum redistribution, star formation (SF) 
and the energy and momentum feedback from the stel- 
lar evolution onto the interstellar medium (ISM). Specif- 
ically, no concensus exists at present on formation of 
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galactic disks, bulges, bars and supermassive black holes 
(SBHs, e.g.. White & Rees 1978; Fall & Efstathiou 1980; 
White 1996; EUis 2000). Because these structures are 
made from the baryonic matter, either gas or stars or 
their mixture, a certain degree of dissipation must be in- 
volved. Detection of quasars at redshifts z ^ 6 (e.g.. Fan 
et al. 2001; Mathur, Wilkes & Gosh 2002; Bertoldi et 
al. 2003), when the age of the Universe was less than 
a gigayear, rises intriguing questions about the origin of 
their host galaxies within the framework of the hierarchi- 
cal clustering, and is important for the fate of the SBH 
— bulge correlation, i.e., M, — Mbuigc, and the masses 
of the quasar host galaxies, at these z (e.g., Trenti & 
Stiavelli 2007). In this work we investigate how various 
dissipative processes, including the SF and the feedback 
from stellar evolution onto the ISM, affect the galactic 
disk and bulge formation and evolution. For this purpose 
we follow a monolithic collapse of an initial perturbation 
in the DM and baryons within the Hubble flow. 

When exactly did the disks and other structural el- 
ements of disk galaxies, like bulges, bars and central 
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SBHs, assemble? How much evolution has occured since? 
Both orbital and groimd-based telescopes have began to 
tighten the constraints on these issues, but essential de- 
tails still missing and controversies are abound. First, 
within the hierarchical merging scenario, the so-called 
'angular momentum catastrophe' — the angular momen- 
tum (J) distribution of the collapsing baryonic matter 
does not match that of the observed disks in nearby 
Universe (e.g., Sommer-Larsen, Gotz & Portinari 2002; 
Mailer & Dekel 2002; van den Bosch 2002b; Burkert 
& D'Onghia 2004). This distribution is characterized 
by low and high J tails which result both in the over- 
extended disks and very massive bulges. In addition, 
gravitational torques between the collapsing baryons and 
the DM halos lead to disks with radial scalelengths sub- 
stantially smaller than those observed. Under the as- 
sumption of ,/ conservation during the collapse, the spe- 
cific J generated by the random work in the expanding 
Universe is too small to explain the formation of bulgeless 
disks in the late type galaxies (van den Bosch, Burkert 
& Swaters 2001; van den Bosch 2002a). 

The second issue is that of the 'core catastrophe' (e.g., 
Moore et al. 1999) and relates to the DM central den- 
sity distribution. The numerical simulations of DM halo 
formation have shown that these attain a universal den- 
sity profile which can be nicely fitted by a two-parameter 
family (Navarro et al. 1996, hereafter NFW). However, 
observations of disk galaxies and galaxy clusters indicate 
that the DM halos have flat density cores rather than 
NFW cusps (e.g., Flores & Primack 1994; Moore 1994; 
Burkert 1995; Kravtsov et al. 1998; Salucci & Burkert 
2000; Boriello & Salucci 2001; de Blok & Bosma 2002; 
Sand et al. 2004) . The standard picture of slow adiabatic 
collapse of a smooth baryonic component in the centers 
of such halos can only lead to a larger central concentra- 
tion, making matters worse. El-Zant, Shlosman & Hoff- 
man (2001) and El-Zant et al. (2004) have shown that 
a flattening of the central density cusps can be a natu- 
ral outcome of energy deposition by the accreted clumpy 
baryons via a dynamical friction both in galaxies and 
in clusters of galaxies. Alternatively, Weinberg & Katz 
(2002) and HoUey-Bockelmann, Weinberg, Katz (2005) 
have claimed that the cusps can be leveled through the 
action of the stellar bars, but this was strongly disputed 
by Sellwood (2003) and McMillan & Dehnen (2005). 

Third, Abraham et al. (1999) and van den Bergh et 
al. (2002) have claimed that the fraction of stellar bars 
sharply decreases above the redshifts of z ~ 0.5, while 
the analysis of GEMS survey (Rix et al. 2004) of 1600 
galaxies has shown that both the fraction of bars and 
their size and axial ratio distributions remain largely un- 
changed up to 2; ~ 1 at least (Jogee et al. 2004). This 
latter evidence is supported by Elmegreen et al. (2004) 
and Sheth et al. (2003) based on much smaller samples 
of 186 and 4 galaxies respectively. 

Finally, the DM halos resulting from cosmological sim- 
ulations appear to be strongly triaxial, both in their pro- 
lateness and flatness,^ €p and fp ^ 0.5 (e.g., AUgood 
et al. 2006). On the other hand, halo shapes inferred 
from observations in the local Universe are axisymmetric, 

* We define the density prolateness as = 1 — b/a, where b/a is 
the intermediate-to-major axis ratio and the flatness as fp = 1—c/a 
with c/a — the minor-to-major axis ratio in the DM halo 



although individual objects may exhibit some mild de- 
gree of prolateness (e.g., Rix & Zaritsky 1995; Merrifield 
2002). Theoretical arguments supplemented by numer- 
ical simulations have shown that bars are incompatible 
with prolate halos (El-Zant & Shlosman 2002; Berentzen, 
Shlosman & Jogee 2006). When taken in tandem with 
the plethora of stellar bars observed locally (e.g., Mari- 
nova & Jogee 2007; Menendez-Delmestre et al. 2007; 
Knapen, Shlosman & Peletier 2000) and at the inter- 
mediate redshifts of ~ 1, corresponding to the lookback 
time of 8 Gyr (Jogee et al. 2004), these results indicate 
strongly that the halo shapes have evolved only mildly 
during this time period. Furthermore, Dubinski (1994) 
and Kazantzidis et al. (2004) have shown that the halos 
lose some of their triaxiality when baryons are added. In 
the next step, Berentzen & Shlosman (2006) have simu- 
lated the disk growth in an assembling DM halo within 
the 4/1^^ Mpc comoving box. The initial conditions have 
been produced by means of the constrained realizations 
algorithm of Hoffman & Ribak (1991). The halo prolate- 
ness was found to be washed out within the disk radius, 
and the halo flatness to decrease substantially. 

The strongly nonlinear process of galaxy formation 
coupled with a long list of dissipative processes provide a 
strong insentive for a numerical approach. A number of 
numerical methods to follow up the disk formation have 
been attempted. In the models of Immeli et al. (2004) , 
when the gas in the proto-disk cools efficiently, the disk 
fragments and forms a number of massive clumps of stars 
and gas. These clumps spiral in to the galactic center and 
merge there to produce a strong starburst, as also shown 
by Shlosman & Noguchi (1993), leaving a central ob- 
ject which resembles a classical bulge. The disk growth 
was simulated in a fixed spherically-symmetric external 
potential of the DM halo by adding the gas directly to 
the disk. Alternatively, Samland fc Gerhard (2003) have 
modeled a disk and the associated DM halo growth by 
adopting the mass inflow rate from separate cosmologi- 
cal simulations, while preserving the spherical symmetry 
of the halo, using a chemo-dynamical code. A number 
of issues have been investigated, among them the stellar 
population of the bulge. A possibility of forming bars 
from baryonic spheroids immersed in triaxial DM halos 
was considered by Gadotti & de Souza (2003; see also 
El-Zant et al. 2003). 

Sommer-Larsen, Gotz & Portinari (2003) analyzed the 
formation history of two galaxies with circular velocities 
comparable to the Milky Way (MW) and found gas accre- 
tion rates, and hence X-ray halo luminosities, at z ^ 1, 
6-7 times larger than at z = for these objects. More 
generally, it is found that the gas infall rates onto these 
disks are nearly exponentially declining with time, both 
for the total disk and for the "solar cylinder." This result 
supports the exponentially declining gas infall approxi- 
mation often used in the chemical evolution models. 

Okamoto et al. (2005) has investigated the effect of the 
SF feedback on the morphological evolution of galaxies 
by assuming that the SF can proceed in two modes — 
either in the high density or shock-triggered fashion. The 
former mode leads to elliptical galaxies while the latter to 
disk ones. The authors made a number of assumptions 
that require a further validation, e.g., the SPH kernel 
was fixed at the constant mass and not on the constant 
number of particles. 
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The simulations of Govcrnato et al. (2004, 2007) form 
rotationally-supportcd disks with rcaUstic exponential 
scalelengths and fit both the / band and baryonic TuUy- 
Fisher relations (see also Stinson et al. 2006; Kaufmann 
et al. 2007). An extended stellar disk forms inside the 
MW-sized halo immediately after the first merger. The 
combination of the UV background and supernova (SN) 
feedback drastically reduces the number of visible satel- 
lites orbiting within this halo, bringing it in a fair agree- 
ment with observations. Formation of stellar populations 
in the MW-typc galaxy were studied by Scannapieco et 
al. (2006) using a combination of high- resolution iV-body 
simulations with semi-analytical methods of metal en- 
richment. Substantial constraints on the formation of 
Population III stars have been obtained. Brook et al. 
(2007b) have followed up with chemodynamical simula- 
tions, including the SF, concluding that surveys of low- 
metallicity stars in the Galactic halo can directly con- 
strain the properties of Population III stars. Further- 
more, the disk evolution was simulated by Brook et al. 
(2007a), examining the possibility that high a elements- 
to-iron abundance ratios detected in the thick-disk com- 
ponent of the Milky Way can be explained by SF during 
the gas-rich mergers. 

All the above studies have brought a very useful in- 
sight into the formation of disk galaxies and have con- 
tributed to a better imdcrstanding of specific aspects of 
the problem. Yet, as already mentioned, the problem 
is very complex and much more remains to be investi- 
gated, in particular regarding the dynamics of forming 
disks, the interaction between the various morphological 
components and subcomponents, as well as their origin 
and properties. The present paper aims to advance of 
some of these open questions. Our emphasis is more on 
the dynamical aspects and includes discussions of the for- 
mation of bars, both inner and outer, made in such ah 
initio simulations. Our algorithm includes the SF, the 
feedback from the stellar evolution, etc. We also per- 
formed a relatively large number of simulations, so as to 
gain insight on how the SF and energy feedback affect the 
structural properties of the galaxy components and sub- 
components. We use cosmological initial conditions and 
follow the Hubble expansion and the subsequent collapse 
of a single perturbation in the gas and the DM in order 
to analyze the role of different processes in the formation 
of disk and bulge components. We aim at understanding 
the role of the external cosmological and internal dissi- 
pative factors. 

Section 2 describes the details of oiu' method of mod- 
eling the SF, our approach in determinig the shapes of 
spheroidal and disk components, and the initial condi- 
tions used for numerical simulations. Sections 3 and 4 
present results of our simulations and model comparison. 
Discussion follows in the last section. 

2. NUMERICAL METHODS 

2.1. Modeling the Gas, Stars and DM 

The numerical simulations have been performed with 
an updated version of the FTM-4.4 hybrid c;ode which 
evolves the coUisionless (iV-body) and gaseous (Smooth 
Particle Hydrodynamics, SPH) components (Heller & 
Shlosman 1994; Heller 1995), with iV ~ 5 x 10^ and 
A^SPH ~ 5 X 10^. In addition, a number of models 
have been evolved using the parallel version, FTM-4.5, of 
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the code. The gravitational forces were computed using 
the public version of the routine f alcON (Dehncn 2002), 
which is about ten times faster than optimally coded 
Barnes & Hut (1986) tree code. The tolerance parame- 
ter is fixed at 0.55. 

The energy equation describing the internal state of 
the SPH particles have been used to evolve a multiphase 
ISM by calculating its heating and cooling rates and it- 
erating for the temperature. Star formation was used to 
convert SPH particles into coUisionless particles which 
exert an energy feedback onto the surrounding gas. The 
gravitational softening for the coUisionless DM and stel- 
lar particles is e*,grav = 150 pc, except for two mod- 
els specified in Table 1. For the SPH particles we use 
two options to soften the gravity — constant soften- 
ing of Eg.grav = 150 pc Or 250 pc, and a dynamic soft- 
ening (Heller & Shlosman 1994) with the minimum of 

^g,min,grav — 250 pC. 

A large number of tests has been performed to check 
the sensitivity of the results to the algorithm and its 
parameters. A typical run with the DM only conserved 
the energy to within 1% and J to within 0.1%. 

2.2. Star Formation 

We anticipate that the SF processes and the associated 
feedback from stellar evolution coupled with an overall 
dissipative and coUisionless dynamics within the forming 
galaxy determine its evolutionary path. Strictly speak- 
ing, because the smallest spatial scales resolved here are 
substantially larger than scales which are relevant for an 
'individual' SF event, the imdcrlying collapse should be 
treated witin the context of the two-fluid Jeans instabil- 
ity (Zel'dovich & Novikov 1975) in a sheared medium. 
We adopt the modified prescription for the SF from 
Heller & Shlosman (1994) that required the gas to be 
Jeans unstable and subject to the feedback from OB stel- 
lar winds and supernovac (SN), and introduced a number 
of statistical elements. Our main assumption is that the 
SF takes place in regions which are contracting, i.e., 

V-v<0, (1) 

and are Jeans unstable, i.e., 

TcoU < ''sound) (2) 

which assures that the collapse would continue unhin- 
dered down to the smallest resolved spatial scales. Under 
this condition, the pressure gradients (which are estab- 
lished over the sound crossing time Tgound) will not be 
able to build up during the collapse time TcoU. Moreover, 
if the gas cooling time satisfies TcooI ^ icoU, this will ease 
the condition given by eq. (2), which, for simplicity, can 
be replaced by pgas > Pcrit = 7x10^^^ g cm^'^ for temper- 
atures in excess of 10^ K (e.g., Navarro & White 1993). 
This happens because the cooling function A(T)/T has 
a minimum at T ^ 1.6 x 10® K. 

For larger densities than pcrit, the short cooling time 
Tcooi will be guaranteed if the neutral (cold) gas will be 
at least moderately self-gravitating in the background 
density pback, 

> QicritPback- (3) 

Here pngas corresponds for the neutral gas (i.e., HI and 
He I), and pback = Phgas + Pstar + Pdm combines the mass 
of the ionized (hot) gas, stars and DM, with acrit ^ 1 
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being the fudge factor. For the standard model, acrft 
0.5 (sec Tabic 1). 

Gas is converted into stars at a rate 



-'ngas 



pngas 
TcoU ' 



(4) 



where the SF (collapse) timescale is taken as Tcoh = 

Q!ffTdyn, with Tdyn = -^/Stt/ 1 6pngas • In reality, the col- 
lapse time for molecular clouds is extended beyond the 
free-fall time, e.g., because of an additional support from 
the MHD turbulence. We have used a.a = 10 for the 
standard model. 

Next, we introduce the probability that a gas particle 
of mass nig produces a stellar particle of mass rUs during 
a given timestep of length At of 



PSF = V 



1 — exp 



fngasAt 
TcoU 



(5) 



where /„gas = Pngas/Pgas and T] = msp/mg. Here 
msF = J^g (1 — Rsf) represents the stellar mass which 
is produced in the SF event. However, a stellar particle 
with a mass of only nis = 0.25mg is created as i?sF = 0.4 
of the mass is instantaneously recycled back to the parent 
gas particle, along with an increment in its metallicity by 
Azg ~ ZsFfnsF, where Zg-p = 0.00793 is the metals yield 
(Leitherer, Robert & Drissen 1992). The quantity 77 is 
not allowed to be less than unity, in which case we set 
T] = 1 and nis = rUg, and the gas particle is removed from 
the simulation. Each gas particle is capable of creating 
four generations of stars with different metallicities and 
Trig decreases accordingly after each SF event. The mass 
of the stellar particle represents a cluster of stars. The 
fraction and timing of massive stars that produce the 
OB line-driven winds and SN are calculated assuming 
the Salpeter IMF. 

2.3. Feedback from Stellar Evolution 

The thermal balance in the interstellar gas is calculated 
using the energy equation. Sources of heating and cool- 
ing include adiabatic, viscous and radiative processes. 
The balance of the specific internal energy along with the 
gas ionization fractions (H, He, He+, He-I— e) and 
the mean molecular weight arc computed as a function 
of density and temperature for an assumed optically thin 
primordial composition gas (Katz et al. 1996). This al- 
lows for the computation of the neutral gas density Pngas 
required for the SF criteria of eqs. (1 3). 

To model the energy feedback from the stellar evolu- 
tion onto the gas, the newly formed stellar particles inject 
energy from SN and OB stellar winds into the A^sf = 16 
surrounding gas particles. The total rate of energy de- 
position by stars, Esf = i^wnd + -E^sn, is given by, 

(6) 



Ey/nd = espPwudmsFm^ z, 



for ts <6x 10^ yrs, and 

-E^SN = eSF-PsNTOSF (7) 

for 3x 10^ yrs < is < 3 x 10^ yrs. Here esF is an energy 
thermalization parameter, tg is the age of a stellar parti- 
cle and Pwnd = -PsN = 2.75 x 10^^ erg yr-^ are the 
energy deposition rates per solar mass by massive stars 
and SN, respectively. The metallicity of a stellar particle 
is given by in units of the solar metallicity (Maeder 



1992, 1993). The timestep of such 'active' stellar parti- 
cles (and of all the gas particles) is restricted to being 
smaller than Atfh = 6x 10^ yrs in order to properly re- 
solve the feedback timescale of Tfb = 3 x lO'' yrs. 

The NsF gas particles in the vicinity of the 'active' 
stellar particles are treated in the following way. Their 
radiative cooling is temporarily disabled when receiving 
the energy from a stellar particle and when the condition 
in eq. (2) is still fuUfilled. The energy ii^sF is thermalized 
and deposited in the gas in the form of a thermal energy, 
then converted to kinetic energy through the equations 
of motion. This method is preferable over injecting a 
fraction of the stellar energy directly in the form of a 
kinetic energy this approach is ambiguous and the 
results are sensitive to the value of a timestep and to the 
particle number on the receiving end. 

For a comparison, we have also run some models with 
an isothermal equation of state (EOS), implementing the 
feedback by directly over-pressuring the gas through an 
increase of the local sound speed over the isothermal 
value. Specifically the sound speed, Cg, is increased by. 



Acg = ^sF^fb (7-1) /2cg. 



(8) 



2.4. Halo Shape Determination 

To determine the intrinsic shape of the DM and stellar 

particle distributions, we remove any residual net veloc- 
ity from the system, reject any unbound particles and 
iterate this procedure to locate the density center (e.g., 
Aguilar & Merritt 1990), 



Td 



(9) 



where pi is the local mass density at the position of each 
particle. 

In the reference frame defined by the so determined 
velocity centroid and density center, we compute the 
eigenvalues of the moment of inertia tensor for the mass 
within a specified radius. From these we may determine 
the axes a > 6 > c of a uniform spheroid with the same 
eigenvalues. The ratios of these axes may then be used 
to characterize the shape of the system, as defined by 



h 
a 

c 
a 



61 — 62-1- 63 



ei + 62 - 63 



A 



(10) 



(11) 



where A = 62 — 61 -I- 63 for the eigenvalues 63 > 62 > 61. 

2.5. Initial Conditions 

The initial conditions are those of a spherically- 
symmetric density enhancement in the Einstein-de Sitter 
Universe. We use an open CDM (OCDM) model with 
Qq = 0.3, h = 0.7, where Oq is the current cosmological 
density parameter and h is the Hubble constant normal- 
ized by 100 km s'^ Mpc"^ Models with fio = 0.3 and 
A ~ 0.7 in ACDM Universe have been run as well, but 
will be discussed elsewhere. Because the collapse time 
of our DM halos happens at 2; ~ 2 (Table 1), the dif- 
ferences between the OCDM and ACDM models are in- 
significant. We have also rerun our standard model N3 
with the WMAP3 baryon fraction of 17% (model N41) 
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TABLE 1 
MODELS 



Moaei 


Zc 


<:SF 


as 


"crit 


^* ,grav 


Eg.grav 




IN otes 


N3 


2 


0.3 


10 


0.5 


0.015 


0.015 




standard model: DM halo+baryons 


QN3 


2 


_ 


_ 


_ 


0.015 


_ 




as N3 but DM only 


N4 


2 


0.3 


10 


0.5 


0.015 


0.015 




enforced axisymmetric halo 


N 5 


2 


0.3 


10 


0.5 


0.015 


0.015 


isothcrma.! 




N6 


2 


0.05 


10 


0.25 


0.015 


0.015 






N 8 


2 


0.1 


10 


0.25 


0.015 


0.015 






N 14 


2 


0.3 


1 


0.5 


0.015 


0.015 






N 15 


2 


0.1 


1 


0.5 


0.015 


0.015 






N16 


2 


0.1 


1 


0.25 


0.015 


0.015 






N17 


2 


0.1 


10 


0.5 


0.015 


0.015 






N18 


2 


0.3 


10 


0.5 


0.025 


0.025 






N19 


2 


0.1 


1 


0.5 


0.025 


0.025 






N20 


2 


0.1 


1 


0.5 


0.015 


0.025 




as N19, with Afsp^noighb = 32 


N21 


2 


0.1 


1 


0.5 


0.015 


0.025 




as N19, with dynamic softening 


















in the gas 


N22 


2 


0.1 


1 


0.5 


0.015 


0.025 


isothermal 


as N21 


N 23 


2 


0.01 


1 


0.5 


0.015 


0.025 




as N21, with esF = 0.01 


N24 


2 


0.01 


1 


0.5 


0.015 


0.025 


isothermal 


as N22, with esp = 0.01 


N25 


2 


0.05 


1 


0.5 


0.015 


0.025 




as N21, with esF = 0.05 


N26 


2 


0.05 


1 


0.5 


0.015 


0.025 


isothermal 


as N25 


N27 


2 


0.3 


1 


0.05 


0.015 


0.015 




as N14, with Ocrit = 0.05 


J\ 28 


2 


0.3 


10 


0.5 


0.015 


0.050 




as N3 


N41 


2 


0.3 


10 


0.5 


0.015 


0.015 




as N3, with WMAP3 baryon 


















fraction 


N42 


2 


0.3 


10 


0.5 


0.015 


0.015 




as N3, with 5 times more 



baryon and DM particles 



Note. — Columns: (1) model (see text); (2) estimated collapse redhift for the DM halo; (3) the fraction of thermalized energy from the 
stellar feedback; (4) the collapse time for SF clouds in units of the free-fall time; (5) fudge factor for the critical density when triggering 
the SF; (6) stellar gravitational softening; (7) limiting value of dynamic gravitational softening in the gas; (8) notes 



and (separately) with 5 times larger number of the SPH 
and collisionless particles (model N42) — the differences 
with the original model appear to be only quantitative. 

The initial density profile corresponds to the average 
density around a 2a peak in a Gaussian random density 
field with power spectrum P{k) = Ak~^ exp (^—k^Rf) . 
The simulations start at Zj = 36 to ensure they begin 
in the linear regime. The particles are initially moving 
outward with the Hubble flow, with consecutive shells of 
mass stopping, turning around and falling back inward. 
The model is parameterized by the filter mass Mf , which 
defines the mass contained within a sphere of radius 2 R{ 
(e.g., Thoul & Weinberg 1995). 

The parameters Mf and Rf are determined by choosing 
the collapse redshift Zc at which the shell Rf reaches the 
center in a non-dissipative collapse and a circular veloc- 
ity Vc associated with the shell surrounding the charac- 
teristic filter mass Mf if virialized at half its turn-around 
radius. In addition, the model is characterized by the 
virial coefficient q = 2Trnd/ where Tmd is the kinetic 
energy in random motions and U is the potential energy, 
and a spin parameter A = J G~^M~^/^, where J, 
E, and M are respectively the total angular momentum, 
energy, and mass of the system. 

The initial angular velocity is taken as lu cc , where 
r is the cylindrical radius and the central kpc has been 
softened. The angular velocity has been normalized by 
the requirement that the initial A = 0.05 for all models. 
This distribution of the specific angular momentum com- 
plies with the universal profile of angular momentum in 
DM halos, based on the statistical sample of halos drawn 



from high-resolution numerical simulations (Bullock et 
al. 2001; see also Fall & Efstathiou 1980). 

For all the models we adopt the values of = 2, g = 
0.01 and A = 0.05. In addition, we have run models 
with Zc = 1 and 6, but those are not shown here. The 
initial models consists of collisionless DM particles with 
the total mass of 7.03 x 10^^ M0 and of the gas, given by 
the SPH particles, which have the same radial profiles as 
the DM. The initial gas comprises 10% of the total mass. 

3. RESULTS: MODEL EVOLUTION 

Our results are presented in this and next sections. We 
describe the evolution starting with QN3, the pure DM 
model, and continue with models which include all com- 
ponents, as well as the SF and the feedback from stellar 
evolution — the standard model N3 is followed by mod- 
els from Table 1. Each siibsection emphasizes a different 
aspect of evolution. We provide a more specific com- 
parison between models, based on parameter variation, 
in Section 4. Note that the specific times characterizing 
the model evolution are given for comparison only. 

3.1. Pure DM Model 

Evolution of pure DM halos within the framework of 

a monolithic collapse is limited to a virialization and in- 
complete violent relaxation in the varying background 
gravitational field. Wc find that the angular momen- 
tum J incorporated in the initial conditions is essentially 
channeled into the internal circulation. In all models, 
including those with baryons, the halo forms a prolate 
triaxial figure and its minor axis is directed along the 
J-axis. The halo figure tumbles with a pattern speed 
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of ilh ^ 1.4 km s^^ kpc^^ during the collapse period of 
t < 1.9 Gyr and ~ 0.2 km s^^ kpc^^ thereafter. This is 
more than two orders of magnitude lower than the pat- 
tern speeds of other prolate figures, i.e., stellar or gaseous 
bars, forming within the growing disk, even after the bar 
secular slowdown. Overall the halo figure makes ^ 200°, 
or less than one full rotation, in a Hubble time. Hence, 
for all practical reasons, we can consider that the halo 
figure is static with respect to any developing bars in the 
system. The initial A = 0.05 parameter remains constant 
for QN3 and shows an increase by ~ 30% for models with 
gas, over the Hubble time. The DM density profile in 
QN3 quickly achieves the NFW shape and the charac- 
teristic radius i?s is ~ 9 kpc by the end of the run. 

Next, we determine the evolution of the DM halo shape 
in QN3. The inner 30 kpc halo turns around first and 
develops a triaxial structure between ~ 0.5 — 1 Gyr. 
Figs, la, c show the evolution of axial ratios, b/a and 
c/a, in QN3. Beyond ~ 50 kpc, the DM collapses over 
^ 5 Gyr and forms a strongly triaxial figure as well. The 
initial triaxiality develops in response to a short-lived ra- 
dial orbit instability (e.g., MacMillan et al. 2006). The 
DM halo in QN3 exhibits a clear increase in its triaxial- 
ity with radius, from the inner, ^ 20 kpc, to the outer, 
^ 40 kpc, halo, both in its prolateness and flatness. The 
former range is 0.1 ^ ep < 0.4, while the latter lies within 
0.35 ^ /p ^ 0.55. Both parameters show a substantial 
increase into the outer halo. 

Roughly speaking, the inner halo is assembled dur- 
ing the first Gyr while the outer one during ~ 5 Gyr. 
Fig. 2 shows the growth of the halo virial mass, Mvir 
which is computed by requiring that the mean halo den- 
sity is equal to the critical density of the universe times 
Ac, where Ac is obtained from the top-hat model (e.g., 
Romano-Diaz et al. 2006, 2007). The halo virial mass 
reaches 7.4 x 10"^^ M© by z = 0. The halo growth rate 
agrees well with the fully cosmological simulations in 
OCDM universe (e.g., Romano-Diaz et al. 2007). Unlike 
Samland & Gerhard (2003), we do not assume the law 
of halo growth and do not limit the shape of the halo to 
be spherically-symmetric. 

3.2. Evolution of the Standard Model and Beyond 

Next, we present models that include baryons and their 
associated processes. The standard model, N3, is treated 
as a representative one, and is compared to other models. 
The A parameter is growing from 0.05 by ^ 30% in all of 
these models. We note that the very slow tumbling of the 
halo figure is not affected by the presence of the baryonic 
matter — the angular momentum is channeled mainly 
into the internal circulation. 

The overall trend in the evolution of a DM halo remains 
the same in the presence of baryons (Fig. la-d). but its 
prolateness and flatness are lessened by ^ 10% — 20% 
compared to that in QN3. This means that the outer 
halo stays substantially more prolate and flattened than 
the inner one, as in QN3. The maximal prolateness and 
flatness of the outer halo are reached earlier in N3 and its 
subsequent reduction is stronger there over the next 2- 
3 Gyr. The inner 10 kpc also become more axisymmetric 
than in QN3. Asymptotically, N3 is substantially less 
prolate than QN3. The washing out of the prolateness 
in the inner 20 kpc is clearly connected with the buildup 
of the disk over the first 2-3 Gyr (as shown already by 
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Fig. 1. — Evolution of the axial ratios b/a (a and c) and c/a (b 
and d) for pure DM halo model, QN3 and for the standard model 
N3. The curves correspond to 2 kpc (solid red), 5 kpc (dashed 
green), 10 kpc (dot-dashed blue), 50 kpc (solid orange), 100 kpc 
(dotted black) and 300 kpc (dash-dot-dot-dotted magenta). 
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Fig. 2. — QN3 model: assembling the DM halo within its 
virial radius. Shown is Mvir/MuM ratio, i.e., the virial mass nor- 
malized by the total DM in the simulation. The inflow rate is 
~ 600 M0 yr"-*^ during the first ~ 1 Gyr and ~ 15 Mq yr~^ aver- 
aged over the next 9 Gyr. 



Berentzen & Shlosman 2006). 

Two factors play the major role in removing the pro- 
latenes of the inner halo. First, the baryonic (density) 
response to the halo (equatorial) prolateness is perpen- 
dicular to the halo's nearly non-rotating major axis. In 
this case the inner inner Lindblad resonance (inner ILR) , 
if it exists, moves to the center, while the outer ILR to the 
very large radii, outside the halo, securing this response 
everywhere within the halo. This baryonic response di- 
lutes the prolate potential of the DM (more about this 
response in section 3.6). Second, the disk sitting deep 
inside the DM potential, as well as numerous baryonic 
clumps forming prior to the disk assembly, scatter the 
DM particle orbits, further reducing its triaxiality. The 



Structure Formation in CDM Halos 



7 



gravitational quadrupole of the disk, during its assembly, 
interacts with that of the prolate halo this induces a 
chaotic behavior in the halo and in the disk. The result 
is a reduced prolateness and flatness in the halo and a 
more axisymmetric disk. 

For the stellar component, the inner ~ 5 kpc of N3 
show a mild prolateness asymptotically, with the axial 
ratio of 6/a ^ 0.95 at ^ 10 kpc due to a mild bar. The 
ability of the DM to maintain the non-negligible prolate- 
ness affects both the dynamical evolution of the galactic 
disk and facilitates the angular momentum redistribution 
in the disk-halo system, as we discuss in section 3.8. 

Because dissipative processes dominate the formation 
of a stellar disk, it is dramatically flatter than the DM. 
For N3 and other models with baryons discussed here, 
the disk formation and evolution can be divided roughly 
into two stages: gas-dominated and star-dominated. The 
former stage lasts over the first 2 3 Gyr — we observe 
that the gas dominates over the stars even in the central 
regions of the disk. The latter stage lasts over the rest of 
the simulation time, during which the ratio of baryonic- 
to-DM matter within the central 10 kpc reaches ^ 0.6 — 
0.65 and subsequently rises to ~ 0.74. The inner few kpc, 
however, are baryon dominated. The gas rotation curve 
flattens at about 290 km s~^. The vertical thickness of 
the gas layer is large over the first 5 Gyrs and then is 
gradually reduced. 

The initially gas-dominated disks are not axisymmet- 
ric, and exhibit a range of morphologies: the m = 2 
and m = 3 modes dominate them at various times, and 
occasionally they are affected by the off-center, m = 1 
mode. With time, the disks become more axisymmetric, 
but some stay oval for a long period of time. 

Even towards the end of the simulation in N3, only 
a fraction of baryons resides within the central 10 kpc 
(where most of the disk lies): 5.8 x 10^" Mq in stars and 
5 X lO'' M0 in gas. The total mass of baryons within this 
radius is 6.3 x lO^*' M©, or 49% of the baryons in the 
simulation. A similar baryon fraction, ~ 40% — 60%, is 
characteristic for all models. Between a third and a half 
of all baryons reside outside the SF region in the disk 
plane (beyond 15 kpc) and in the DM halo, i.e., in the 
form of partially rotationally-supported and/or hot gas. 
The rest of the gas is shock-heated to virial temperatures 
of ~ 10^ K, and ionized. 

The DM density profile for N3 can be fitted by the 
NFW law although the quality of the fit is lower than 
that of QN3, and the final value of Rg is ^ 4.4 kpc. The 
inner few kpc and the DM cusp are dominated by the 
baryons at the end of all the runs with the gas. 

3.3. Early Baryon Inflow and the Disk Formation: 
the Cat's Cradle 

Early baryon inflow towards the center leads to the 
formation of a speciflc configuration near or at the DM 
density peak that resembles the "cat's cradle" pulled by 
a number of nearly radial 'strings' (filaments). This is 
observed in all the models, except in N4 which has an 
enforced axisymmetric halo. The source of this mor- 
phology is the background DM distribution within the 
inner ~ 30 kpc halo during the initial ~ 1 Gyr (Fig. 3, 
left). The baryonic strings have one-to-one correspon- 
dence with the DM strings. Furthermore, the face-on 
disks show such strings attached to their edge, where the 



main baryonic inflow joins the disk. The strings arc not 
perfectly radial, and at the disk edge arc tangent to it — 
the cold gas inflow, therefore, joins the disk without be- 
ing shocked, thus supplying a cold, unvirialized gas to the 
disk. The typical symmetry resembles the irregular poly- 
gons, and ranges from triangles to octagons (e.g.. Fig. 3, 
right). The central object — an amorphous clumpy pro- 
todisk, has a typical size of 2-4 kpc and the radial strings 
attached. The protodisk has an elongated shape and its 
figure initially nearly stagnates. The baryonic blobs tend 
to form on these strings in the equatorial plane and fall 
along them onto the disk, visibly speeding up its figure 
rotation. The feeding of the disk along the strings con- 
tinues for some time during the first Gyr. Models N23- 
N26 with a small SF feedback, esF = 0.01, exhibit larger 
clumps. The early disk appears as a fat configuration 
with an equatorial SF and with flaring edges. By about 
~ 1 Gyr, the disk grows substantially, its shape changes 
rapidly between a less and a more prolate one, and with 
more or less prominent spiral structure. Dynamical fric- 
tion appears to play a major role following the clump 
joining the disk similar to the behavior of the gas-rich 
disks analyzed by Shlosman & Noguchi (1993). 

3.4. Later Stage: Baryon Mass Inflow History and the 
Disk Growth 

With the beginning of the collapse, the rate of bary- 
onic inflow into the inner 20 kpc reaches its maximal 
rate within the first 1.5 Gyr and then steadily falls off. 
The SF rate shows a similar behavior (except in N27), 
with a peak at 1.5-1.8 Gyr of ~ 38 Mq yr~^ for the 
N3 disk (Fig. 4). Overall, the range of the SF rates is 
10 — 60 M0 yr^^, and the peak SF rate can be delayed 
by 1-2 Gyr compared to the standard model. N27 differs 
and shows a flat SF rate for the first 6 — 7 Gyr. 

The SF rates decay exponentially, in nearly all mod- 
els (except N27), to ^ 5 Mq yr^-"^, at ^ 5 Gyr, and 
^ 1 Mq yr~^ after ~ 10 Gyr. The latter rate is typi- 
cally observed in the normal disk galaxies of the nearby 
Universe. 

In all models, the proto-disk plane has been established 
and visible already at r ^ 0.5 Gyr, across the inner 
12 kpc, while the stars exist only at the very center by 
that time. Thereafter, a small r 2 kpc gas disk be- 
comes prominent and doubles in size to 5 kpc in another 
0.5 Gyr, with the SF at its midplane. After ~ 1 Gyr, 
a gaseous bar of r ^ 2 kpc forms, supplemented by an 
ongoing SF, but collapses to the center and is damped 
shortly (more about these bars in sections 3.5-3.7). By 
T ~ 1.2 Gyr the disks grow almost to 10 kpc with an 
intense SF — the gas component remains puffed up sig- 
nificantly. 

By r ~ 2 Gyr, the gas layer starts to cool down. 
The transient spiral structure becomes visible and is de- 
lineated by the SF regions. By ~ 6 Gyr, the vertical 
structure of the stellar disks in most of the models, but 
not in all, looks alike. We observe that the disk cools 
down and thins substantially when the SF rate drops be- 
low ^ 5 Mq yr^^ this is true for all the models, and 
typically happens after about 5 Gyr. 

All the disks, except that of N27, appear to form from 
inside out, and are gas-dominated initially. The typical 
time for the disk buildup, i.e., when it reached about 50% 
of its final mass at 2: = 0, is 2 Gyr. The stellar rotation 




Fig. 3. — Snapshots of inner DM halo (loft) and disk (right) formation in the N3 model. The time, t = 0.71 Gyr was selected to 
emphasize the 'cat's cradle' morphology of the proto-disk which displays the attached radial filaments. This morphology is driven by an 
underlying DM potential of the forming inner halo. The disk particles consist of gas (green), SF particles (blue) and stars (yellow). Left 
frame size is 30 X 30 kpc and right frame is 12 X 12 kpc 
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Fig. 4. — Star formation rates in the N3 and N27 model disks. 



curve gradually extends to larger radii and steadily rises 
with time, while remaining flat in r. Its peak moves in 
to smaller radii, then recedes somewhat (Fig. 5). (See 
also the comment on the bar pattern speeds evolution, 
reflecting the growing mass concentration in the disk, in 
section 3.6.). The gaseous component, in the inner few 
kpc, forms a high surface density 'thin' disk embedded in 
the thicker 'corona.' This corona extends to larger radii. 
The extent of the stellar disk is limited to the narrow gas 
disk only. While the relative details of this disk-corona 
morphology differ from model to model, the overall struc- 
ture is robust. Models N19 to N26, all show larger disks, 
some with a strikingly visible bulge accompanied by SF 
activity in the nuclear rings. 

The central kpc of an assembling galaxy in N3 is dom- 
inated by the DM during the first ^ 1 Gyr, and subse- 
quently by the gas (for a short period of time) — then by 
the newly formed stars. The baryon/DM ratio tends to 



0.8 within the disk radius at the end of the simulation. 
This value, however, is hardly a standard. It is affected 
by the feedback from the stellar evolution, by the critical 
density for SF and by the collapse time of unstable SF 
clouds. 

3.5. Early Nuclear Bars 

Early gas-dominated (or gas-rich) nuclear bars appear 
as a common feature after about 1 Gyr in all of the mod- 
els, albeit their longevity differs from model to model 
(typically 4-5 Gyr), e.g., N3 in Figs. 6-7. Their strength, 
measured by their axial ratio or by the amplitude of their 
m = 2 mode, changes with time, and they become more 
star-dominated as the gas is channeled inward and pro- 
gressively concentrates in the central ~ 100 pc. The 
(radial) sizes of these nuclear bars change with time as 
well, from few 100 pc to about 2 kpc. The central SF 
is clearly concentrated in these bars. The stellar bulges 
are observed to dominate these bars after 2.5-3 Gyr in 
some models. The vertical extent of the gas layers within 
the bulges appears to be similar to elsewhere in the disk. 
We note that nuclear gas bars, devoid of the stellar com- 
ponent, are observed in the central kpc of our models 
as well. These frequent gaseous nuclear bars survive for 
a limited time, depending on the gravitational coupling 
from the background potential, and collapse to the cen- 
ter. 

Nuclear rings appear prominently around the nuclear 
bars from time to time. They seem to be connected to the 
large-scale disk by pair of gas spiral arms which are also 
visible in the SF (the outer disks are visibly oval), e.g., 
N19 in Fig. 8. These rings become more prominent when 
the gas in the oval disk is channeled to the central kpc. 
The gas nuclear bars, seen to develop inside the ring, 
shrink catastrophically and speed up, until the gas is 
dumped at the very center — leaving a profound central 
hole in the cold gas distribution. 

3.6. Large-Scale Stellar Bars 

Large-scale bars, ~ 3 — 6 kpc, have always developed 
very early, by r ~ 1 Gyr, e.g., in N5 — in all cases they 
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Fig. 5. — Snapshots of a circular velocities in the assembling disk/halo system of model N3 at 1 Gyr (left), 2 Gyr (center) and 12 Gyr 
(right). Total velocity is given by the solid line, DM by dotted line, stellar by dash-dotted, and gaseous by dashed lane. 
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Fig. 6. — Left frame: The evolution of m = 2 Fourier component A2 in the N3 disk. Shown axe stars (left frames) and gas (right) within 
the central 5 kpc (upper frames) and 1 kpc (lower). Right frame: The pattern speed of a stellar bar within the inner kpc of N3 smoothed 
by fitting a cubic spline to the P.A. of the m = 2 in Fig. 60 and differentiated. 



have decayed in a few Gyr. In addition, bars have formed 
late in the evolution, after 5 Gyr, e.g., in N16 (Fig. 9), 
in about 2/3 of the models. The early bars appear in 
response to the prolateness in the DM halo (and of the 
stellar disk consequently) — they do not develop as a 
result of a classical bar instability. These bars are mod- 
erately strong, with e 0.4 — 0.5, and seem to involve 
a large fraction of the stellar disk mass. Their shapes 
as well as the shapes of the outer disks depend on the 
mutual orientation between the DM halo, disk and bar 
major axes. We use the m = 2 mode amplitude, A2, to 
quantify the strength of this mode in the disk and the 
bar, as a function of their P.A. with the halo. Applying 
this procedure to N3, during the first Gyr, shows that the 
disk A2 integrated between 8 kpc - 12 kpc displays its 
minima/maxima when the disk is oriented at — 45°/+45° 
to the halo, i.e., when it leads/trails the DM halo major 
axis. During the second Gyr, the A2 has maxima when 
the disk is oriented at —90° to the halo in about 1/3 of 
the rotations. In other cases, the A2 appear to 'librate' 
around P.A.» 90°. As the disk becomes less oval due to 
its interaction with the halo and the bar, the residual A2 
always librates around P.A.« 90°. 

The same analysis repeated for the bar in N3 shows 
that early in the evolution, the bar 'stagnates' at P.A.~ 
90° to the halo. At later times, the bar sometimes ap- 



pears stronger at this angle. But in most cases, 2/3, 
A2 appears irregular as a function of P.A. Hence the bar, 
as the least massive object in comparison with the disk 
and the halo, shows a much less regular behavior of its 
A2. 

We also find that stellar bars forming late in the simu- 
lation are shorter, than the early bars, e.g., in N16, where 
the bar starts to grow anew after ^6 — 7 Gyr and contin- 
ues its slow growth until the end of the simulations. In 
addition, early bars, large-scale and nuclear, are visibly 
dominated by the SF and are gas-rich. While their late 
large-scale counterparts appear to be gas-poor. 

Fig. 6b exhibits the evolution of the bar pattern speed 
in N3. Due to the disk growth from inside-out, it shows 
rather an unusual behavior — it appears at a fixed, ~ 
90°, angle to the halo and remains in this position for a 
large fraction of a Gyr. Then the bar speeds up sharply, 
till ^ 2 Gyr. By that time the central bar- hosting regions 
in the disk are mostly built and further addition of the 
mass to the disk happens well outside. 

3.7. Nested Bars and Spirals 

In the models presented here, the growing stellar disks 

preserve their strong ellipticity (prolateness) over the 
first 1-2 Gyr. Strong nuclear gaseous bars form and 
experience a distributed SF over their surface. In some 
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Fig. 8. — Example of nuclear rings around nuclear bars: N19 is shown in gas (left frame), gas+stars+SF (right) at r ~ 3 Gyr. 



models, these bars appear tumbling much faster than the 
(oval) disk figure. A pair of grand-design spiral arms ex- 
tend from the surrounding nuclear rings outward, across 
much of the stellar disks. After few rotations, the nu- 
clear gas bar collapses to the center (e.g., N18 and N26 
models in Figs. 10, 11). 

These phenomena become progressively more pro- 
nounced along the sequence N16— >N22, excluding N17 
and N20. N23-N25 are more clumpy initially (corre- 
sponding to Shlosman & Noguchi [1993] scenario with 
clumps spiraling in because of a dynamical friction) and 
with no pronounced nuclear morphology later. The ap- 
pearance of these clumps can be clearly correlated with 
a sharp decrease in the thermalization of stellar feedback 
energy, down to esp = 0.05. 

Model N26 represents one of the clearest cases of nested 
bars having different pattern speeds (Fig. 11). The nu- 
clear bar appears first within the central kpc at around 
1 Gyr as a result of the gas inflow from larger scales, 
where m = 2 and to = 3 are clearly visible and con- 
stantly evolve. This nuclear bar can be observed both in 
gas (and the associated SF) and in stars. At r ~ 1.5 Gyr, 
the nuclear gaseous ring appears with an associated SF 
and is connected to the outer disk by a pair of grand- 
design spiral arms. At this time, the nuclear bar seems 
to have a different pattern speed from these arms and 



tumbles faster. After ~ 3 Gyr, a large-scale stellar bar 
becomes visible within the central 5 kpc. The shape and 
the existence of the nuclear ring depends on the mutual 
alignment of the nested bars (Shlosman & Heller 2002). 
The ring is also visible in the edge-on disk as being dense 
and much thinner than the gas layer within the central 
kpc. we also confirm that the strength of the small bar 
depends on the mutual orientation of the bars — it is 
stronger when the bars are aligned and weaker when they 
are normal to each other (Heller, Shlosman & Englmaier 
2001; Shlosman & Heller 2002; Englmaier & Shlosman 
2004; Debattista & Shen 2006). 

We note that the nuclear bars are gas-rich in the be- 
ginning of this process and form stars vigorously. Grad- 
ually the gas in the bar is dumped on the center, and 
the bar dissolves after r ~ 4 Gyr. The pattern speeds 
of the nuclear bars in all models increase over the first 
1 — 3 Gyr and decline thereafter. The stellar bar dis- 
solves by T ~ 10 Gyr. Nuclear ring survives and becomes 
amorphous. A specific case of evolution of nested bars 
has been analyzed by Heller et al. (2007). 

3.8. Angular Momentum Redistribution 

The pure DM models conserve the total angular mo- 
mentum J with the precision of 0.1%. In models with 
baryons, J is not conserved apriori because of the pro- 
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Fig. 9. — Examples of early (gas-rich) and late (gas-poor) large-scale bars shown in gas-|-stars+SF colors; 8 kpc bar in N5 (left frame) 
a,t T ^ 1.5 Gyr, and 3 kpc bar in N16 (right frame) at the end of the simulation. 
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Fig. 10. — Nuclear bar tumbling with a pattern speed much higher than the oval stellar disk in N18 at r ' 
colors, and shown within the central 12 kpc (left) and 2 kpc (right). 



1.7 Gyr, seen in gas-|-stars+SF 



cesses associated with the SF and feedback. Energy input 
from stars leads to pressure gradients, i.e., hydrodynam- 
ical torques, etc. We find that J of all particles in mod- 
els with SF is nevertheless a well preserved quantity — 
across the models it increases by a small ~ 0.3% — 0.4%. 
These processes, therefore, appear much less important 
than the gravitational torques for the J balance. 

For baryons, J decreases by about 25% — 30% over the 
Hubble time, while for the DM matter J increases by 
about 2.5% — 3%. Most of this decrease happens during 
the first 3-5 Gyr. For stars, J saturates after 5 — 7 Gyr, 
while J in the gas continues to decrease — by the end 
of the simulations about 25% — 50% of baryonic J is in 
the gas, depending on the model. Note that this involves 
the gas at large radii as well, in the regions devoid of 
stars. This explains the large contribution to J by the 
gas compared to its overall low mass fraction at the end. 

Of a special interest is the balance of J within the cen- 
tral 10 kpc, Jio, where most of the baryonic disk resides. 
This angular momentum behaves in the following way 
(e.g.. Fig. 12 for N3). The peak of Jio in the DM is 
achieved at 1.2 — 1.9 Gyr with the follow up sharp 
decrease by a factor of 2-2.5. The central 10 kpc re- 



mains dominated by Jio contributed by the DM for the 
first ^ 3 Gyr. The rise in the baryonic Jio is much slower 
(than for the DM) within this region and it reaches the 
maximum only at ~ 8 Gyr, in some models declining 
about 25% thereafter. The peak of Jio for the DM has a 
similar value, within a factor 2, to that of the baryons 
within 10 kpc (e.g.. Fig. 12). Even more interesting is the 
close similarity between the values of Jio for the baryons 
and DM after the first 'splash ' over most of the evolution 
time — this 'equipartition' of Jio tells that the angular 
momentum transfer from the disk to the halo is very ef- 
ficient and saturates when both components have equal 
J. This behavior is observed for all models. 

Stellar Jio mimics that of the baryons in total, being 
about 25% lower at all times. For the gaseous compo- 
nent, Jio reaches its maximum between 1.5-3 Gyr, within 
a subsequent sharp decline (by a factor of 2) and a slow 
decline thereafter. The gas Jio is always much smaller 
than in the DM (by a factor of 2). The evolution of gas 
J mimics (somewhat) that of the DM overall. 

The specific angular momentum j of the gas is de- 
fined as Jgas/-^^gas, and similarly for stars and the DM. 
In models with baryons we find that jdm increases by 
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Fig. 11. — Differing pattern speeds of nested bars in N26 — two snapshots. Upper: large-scale stellar bar in gas+stars+SF (left), 
nuclear bar in gas+stars+SF (center) and in gas (right) at r ~ 3.00 Gyr. Lower: same as above but at r ~ 3.07 Gyr. Left frame sizes are 
12 X 12 kpc, center and right frames are 2x2 kpc. Note that the bars are nearly perpendicular at the former time and nearly aligned at 
the latter time. 
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Fig. 12. — Evolution of the angular momentum, Jio within 
10 kpc in N3: for the DM (black dotted lane), for the baryons (red 
solid lane), for the stars (blue dash-dotted) and for the gas (green 
dashed). 



about 1% — 2% over the Hubble time, most of this in- 
crease is dated by z 9 — 10. The baryons lose about 
15% — 25%ofj, mainly during the epoch z ~ 4 — 7 . How- 



ever, j of the gas and stellar particles separately increases 
(nearly) always monotonically. Initially, stars form at the 
very center with the minimal j, then the SF propagates 
outward consuming the gas with higher j. Again, the 
specific angular momentum j of the gas increases with 
time in all models because the SF first consumes parti- 
cles with lowest j. Typically, the stellar j saturates its 
growth after few Gyr, while the gas continues to increase 
its j. Note, that the number of gas particles is not de- 
creasing initially, but their mass is — the algorithm is 
using the option of producing generations of stars from 
the same gas particles before removing the gas from the 
balance sheet (e.g.. Fig. 13). 

4. RESULTS: COMPARISON OF MODEL EVOLUTION 

Here we describe the model evolution based on the vari- 
ation of specific parameters related to the SF processes 
and to the feedback from stellar evolution. 

4.1. Varying the Equation of State: Isothermal Models 

The isothermal EOS has been used in N5 and in N22, 
N24 and N26. Compared to N3, all the isothermal mod- 
els develop stronger large-scale (~ 3— 4 kpc size) bars and 
stronger nuclear bars. The large bars drive much more 
pronounced outer spiral arms, both in stars and in gas, 
and accumulate increasingly clumpy gas, especially dur- 
ing the first ^ 1.3 Gyr. (This is also true when compar- 
ing models with and without isothermal EOS: N22 with 



Structure Formation in CDM Halos 



13 



0.5 



Time [Gyr] 



Fig. 13. — Evolution of specific angular momentum j for gas 
(green dashed lane), stars (blue dot-dashed), baryons (red solid) 
and the DM (black dotted) in the N3 model. 



N21, and N24 with N23.) These large bars shrink radi- 
ally (arms persist in the gas) to within the central kpc 
and largely dissolve due to the quadrupole interaction 
with the halos, and to a lesser degree due to the heating 
by a dynamical friction from the numerous clumps, just 
to develop sometime later on. Compared to N3, N21 and 
N23, the bars in the isothermal models show a higher rate 
of SF. The nuclear bars persist for a few Gyrs, leaving 
an oval distortion thereafter. Similarly to the large-scale 
bars, nuclear bars are 'resurrected' from time to time. 

In all of the above models a central opening of about 
1 kpc in the gas appears, surrounded by a ring, by 
r ~ 5 Gyr. Gaseous bars driving two gaseous spirals de- 
velop within this hole, from time to time, and exhibit a 
faster pattern speed than the surrounding flocculent spi- 
ral structure. By t 7.5 Gyr, the hole grows to about 
2-3 kpc in radius and more prominent in the models with 
the isothermal EOS. The comparison of pattern speeds 
shows a diverging evolution for the gas inside 1 kpc and 
5 kpc, e.g., late (after 9 Gyr) speedup of the gaseous 
component in N26, or a late gaseous bar in N22. 

Overall, the isothermal models show morphological 
features that are sharper, including nuclear rings and 
bars, and their SF rates exceed that of the non- 
isothermal. The spiral arms outside the nuclear rings 
have a grand-design appearance. During the maxima of 
the SF, the gas disk is more puffed up. The edge-on 
gas disks exhibit a characteristic concave shape during 
the first 3 Gyrs, and the SF is more pronounced in the 
mid-plane. In comparison, the N3 gas disk is more con- 
vex. The outer massive spirals in the isothermal models 
appear to be responsible for these differences, which dis- 
appers by ~ 3 Gyr. The bulges in the isothermal models 
extends vertically beyond the disk after r ~ 2.3 Gyr, 
and, by the end of the simulations, are somewhat more 
pronounced. The disks in non-isothermal models are 
smaller, e.g., in N3 by about 20%, both in stars and 
gas. Typically, much more gas remains outside the N3 
disk, both in the mid-plane and above the disk. The disk 



shapes differ dramatically during the evolution between 
non-isothermal and isothermal models and the morphol- 
ogy of spiral arms looks different as well. 

4.2. Varying Feedback From Stellar Evolution, esF 

In order to investigate the effect of a stellar feedback 
onto the gas, we have varied esF from 0.3 in N3, to 0.1 
in N17. Additional relevant sequences in Table 1 are the 
'doublets' N8 (esF = 0.1)^ N6 (0.05), N14 (0.3) ^ N15 
(0.1), and the isothermal N22 (0.1) N24 (0.01), and 
the 'triplet' N21 (0.1) ^ N25 (0.05) N23 (0.01). This 
parameter has a major impact on the evolving morphol- 
ogy in our models. The SF starts earlier and tends to 
be more patchy, underlining the clumpy character of the 
gas, when egp is smaller. The clumps spiral in toward 
the center because of the dynamical friction, especially 
around r 1 — 2 Gyr. Models with a smaller feedback 
also show a more prominent morphological features, like 
bars, nuclear rings and spiral arms. 

In the first doublet sequence, the N17 evolution is sim- 
ilar to that in N3, except the SF extends to larger radii. 
The N3 disk is the smallest among the models listed 
above and the N23 disk is the largest. In the early stages, 
the central gas concentration and its central surface den- 
sity are observed to increase with decreasing egp. The 
gas layer thins along the same sequence. 

The stellar bulges form somewhat earlier in models 
with lower esp, in N23 by r 1.6 Gyr, in N25 by 
- 1.9 Gyr, and in N21 by - 2.5 Gyr. Bulges in N17 
and N3 form by ~ 2.5 Gyr. By the end of the sim- 
ulations, the largest bulge is found in models with the 
smallest esp, among the triplet and doublet parameter 
sequences given above, i.e., in N23 (Fig. 14) and N24. 
The smallest bulges are found in N3, N22. Such a bulge 
is also expected in N21. However, as can be seen in 
Fig. 14, the bulge of N25 is smaller than that of N21. 
The reason for this behavior lies in that N21 develops a 
strong late bar, while N25 exhibits only an oval distor- 
tion after the first 5 Gyr. N23 lacks a late bar totally. 
Stellar bars arc known to heat up the central disk re- 
gions and inflate them vertically by increasing their ve- 
locity dispersions, both dynamically and secularly (e.g., 
Athanassoula & Misiriotis 2002; Berentzen et al. 2007). 
Additional heating exerted by the bar in N21 contributes 
to the appearance of a more prominent bulge there. 

We infer that two processes, at least, emerge as being 
responsible for the bulge and stellar bar evolution in our 
simulations — both of these processes are 'guided' by esp, 
the stellar evolution feedback onto the gas. We observe 
that a decrease in egp alleviates the dumpiness in the gas, 
which in turn heats up the stellar component in the disk, 
therefore, delaying or damping the bar instability there. 
This explains why N21 bulge does not 'fit' in the sequence 
shown in Fig. 14. We conclude, therefore, that by the 
end of the simulations there is a clear anti-correlation 
between the bulge prominence and egp. 

4.3. Varying Critical Density for Star Formation, ttcrit 

We look for the effect of the critical gas density for the 
SF by comparing N14 (acrit = 0.5) with N27 (0.05), and 
N15 (0.5) with N16 (0.25). The main difference observed 
along this sequence is the delayed evolution in N27, where 
^ 3 kpc bar develops much later compared to N3 and 
N14, after 8 Gyr only. It can be also seen in the SF 
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Fig. 14. — Stellar bulges in the edge-on disks at the end of the simulations for models N21 (left), N25 (center) and N23 (right). This 
sequence corresponds to a decreasing value of the feedback parameter, egp (see Table 1). 



colors. The dominant SF in the central kpc of this model 
also appears very late, after ~ 2 Gyr, and spreads out in 
r and in z (low surface and volume densities in the gas). 
It picks up, however, after 5 Gyr, and remains much 
higher than in other models. The stellar disk is clearly 
much hotter and puffy in N27. The bulge develops late 
but becomes more pronounced than in N14. 

On the other hand, using a low value adit — 0.05 in 
N27 produces the smallest r ~ 4 — 5 kpc disk with a very 
high surface density, enveloped in a low density stellar 
and gas material (Fig. 15). The SF is concentrated in 
the disk mid-plane. A bar-like stellar and gas response to 
the halo prolateness is evident initially within the central 
^ 3 kpc, but decays after 3 Gyr. This response forms at 
a fixed, nearly normal angle with the DM major axis. 

Hence, using a low density threshold for SF, in the 
presence of a substantial feedback, heats up the stellar 
and gas components in the growing disk. N27 is prob- 
ably the only model where the SF is distributed and is 
triggered at all radii < 15 kpc in the disk — this disk 
does not appear to grow from inside out. It produces a 
substantial bulge but the smallest stellar disk. 

4.4. Varying Timescale for Star Formation, ag 

A number of models have a much shorter timescale for 
SF, ag — 1 instead of 10 (see Table 1) — we compare 
N3 (as = 10) with N14 ag = 1, N8 (10) with N16 (1), 
and N17 (10) with N15 (1). 

Early enough, by r ~ 0.9 Gyr, the SF in N14 is less 
concentrated toward the midplane and the stellar and 
gas disks are visibly thicker than in N3. The SF is also 
more profound in the N14 bars. The gas surrounding 
the disk in N14 is more noticeable, to the extent that 
the gas disk flares beyond 5 kpc. By ~ 1.4 Gyr, a 3 kpc 
disk of higher surface density emerges in N14 — while 
the vertical thickness of the disk is now comparable to 
N3, the radial size of N14 disk is smaller. By ^ 7 Gyr 
the bulge of N3 seems more pronounced, while toward 
the end the N14 bar is more gas-rich than in N3. 

N15 shows earlier and higher rate of SF and has an 
initial stellar, gas (and SF) bar. This bar weakens and 
strengthens occasionally but after ~ 5 Gyr grows mono- 
tonically. For a comparison, N16 develops such a large- 
scale bar by 7 Gyr. The disk in N15 is more puffed up 
than in N17 and the grand-design spiral arms are much 
less visible. Grand-design spiral structure is strong be- 
tween ^ 2.5 — 5 Gyr — apparently when the bar/disk 
are aligned at certain angles to the DM halo major axis. 
The SF, which delineates the arms, follows them, while it 
is more diffuse in N17. The disk is thick and concave at 
this time, much more than in N17. Compact gas clumps 
are present thereafter. By 3.5 Gyr, the SF in N15 has 
subsided more than in N17, the disk seems to be some- 



what smaller and the bulge is more compact as well. It 
also has less gas in the central kpc. The gas layer starts 
to thin first in this model compared to N17. Towards 
the end of the simulations it is clear that the disk in N15 
is overall thicker, and the bulge is standing out with re- 
spect to the face-on disk more than in N17. This bulge 
has an oval/barlike shape. 

Overall, models with smaller as form somewhat 
smaller in size but more puffed up stellar disks. 

4.5. Standard Model with an Imposed Axisymmetric 

Halo 

To enforce the halo axial symmetry during a continu- 
ous infall of material, we have randomized the azimuths 
of 5% DM particles, every timestep. Therefore, no initial 
bars are either triggered or driven by the halo prolate- 
ness. Furthermore, this should lower the ability of the 
halo particles to be in resonance with the nonaxisymmet- 
ric features in the disk. This, in turn, should decrease 
the efficiency of angular momentum flow from the disk 
to the halo and change the stellar bar evolution. 

Indeed, in model N4, the initial collapse proceeds with- 
out the 'Cat's Cradle' phase and no DM filaments are 
present. Consequently, the gas does not form the mas- 
sive clumps. The disk remains nearly axisymmetric at all 
times and no substantial spiral structure is visible in the 
outer part, unlike in in all other models. Some stronger 
oval distortion develops inside 1-2 kpc after ^ 8 Gyr for 
some time and decays thereafter. The gas bars are vis- 
ible from time to time in this region. This evolution is 
in a sharp contrast with N3 which displays an early and 
late bar activity on various spatial scales. The SF and 
spiral morphology in N4 also differ substantialy. 

4.6. Varying Gravitational Softening, Cgiav 

Models N18 and N19 have larger gravitational soften- 
ing in gas and in stars — Cgrav = 250 pc instead of 150 pc 
in other models. We first compare N18 with N3. The 
general evolution of N3 and N18 is similar but it differs 
in some details. First, the former model shows sharper 
features (shocks, spiral arms, SF regions). Second, the 
gas there is somewhat more clumpy and is able to con- 
centrate more towards the center and the disk mid-plane. 
Already during the first Gyr, N3 gas shows a larger sur- 
face density in the central kpc. Third, the bars that form 
in this model are less strong than in N18 (Figs. 6, 16). 
While the N3 nuclear bar is gas poor (efficiently chan- 
neling it to the center), the N18 bar remains gas-rich for 
a few Gyrs. Two factors may reduce the bar strength 
here: larger gas and stellar concentration to the center 
and more clumpy gas that scatters bar stars. This lat- 
ter feature can also prevent the bar orbits to lock in the 
resonance with the DM halo orbits, thus resulting in a 
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Fig. 16. — The evolution of m = 2 Fourier component A2 in 
the N18 disk. Shown are stars (left frames) and gas (right frames) 
within the central 5 kpc (upper frames) and 1 kpc (lower). 



weaker bar. Finally, the stellar bulge in N3 is less promi- 
nent than in N18 but is more centrally concentrated. 

Next, we compare N21 with N19. The stellar soften- 
ing is the same but the gas softening is 150 pc in N21 
instead of 250 pc in N19. As in the previous pair, the 
morphological features are more detailed in N21. This is 
noticeable in the distribution of SF and gas in the strong 
spiral arms — both are more clumpy, but not in stars. 
As before, N21 develops a more massive gas concentra- 
tion in the center, but the difference between the models 
is not as obvious, apparently because the stellar compo- 
nents have identical softening. The nuclear rings appear 
(vertically) thinner in N21 and are much more long-lived. 
Occasionally, they collapse to the center. In the long run, 
the difference emerges also in the stellar distribution as 
well — the bulge is visibly larger in N19 and less cen- 
trally concentrated. The bars in N19 are also more gas 
rich, but periodically this property is inverted between 
the models. On the other hand, the bars in N19 appear 
stronger. Hence, the gravitational softening emerges as 



an important parameter which drives the evolution in 
our models. Especially in gas, this parameter reflects 
the overall smoothness of the gas and can be related to a 
number of intrinsic processes in the ISM, heating/cooling 
and the topology of the magnetic fields — those affect 
the gas dumpiness. 

5. DISCUSSION AND CONCLUSIONS 

We have simulated the formation of a two-component 
galactic disk within an assembling DM halo, including 
the accompanying SF and the stellar feedback onto the 
ISM. For this purpose, we have followed the collapse of an 
isolated cosmological density perturbation with the spin 
parameter A — 0.05, from its linear regime at z = 36 and 
up to z = 0. Our approach has several advantages. We 
have improved the SF criterion, accounting not only for 
the contraction and the critical density of the SF region, 
but also for the background density. Each gaseous parti- 
cle can make up to four generations of stars, with differ- 
ent metallicities, while an additional parameter allows to 
consider a delayed collapse of the molecular clouds due to 
MHD turbulence. We consider the energy feedback both 
from stellar winds and from SN, and we are particularly 
careful in the way this energy is deposited in the gas. 

We have employed a reasonably high number of par- 
ticles, allowing us to resolve a number of important 
features. In particular, our simulations are fully self- 
consistent, in order to avoid artifacts that could be in- 
troduced, e.g., by a rigid halo (e.g., Athanassoula 2002). 
This enabled us to study interesting dynamics linked to 
the halo shape. We have voluntarily refrained from us- 
ing a yet higher number of particles, in order to be able 
to run a large number of simulations, much more than 
is usual in similar studies. With our ^ 30 models, we 
tried to gain insight on how the various parameters, and 
particularly the SF and the energy feedback, affect the 
structural properties of the bulge, the disk and its com- 
ponents, like bars and spirals. Of course, the available 
parameter space is very large and much further work and 
many more simulations are necessary before any thor- 
ough understanding is reached. Yet, even with the num- 
ber of simulations we have at our disposal, several clear 
trends emerge. 

In all models the disk forms over a period of time of 
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2 — 3 Gyr, i.e., it reaches 50% of its final mass over 
this time, although its equatorial plane is established 
very early. The DM halo develops a characteristic tri- 
axial shape, i.e., prolate in the plane perpendicular to 
the original angular momentum, and flattened along the 
rotation axis. For pure DM models, the axial ratios of 
the inner (< 20 kpc) halos lie between b/a ^ 0.8 — 0.9 
and cja ^ 0.75 — 0.85, and of the outer (> 50 kpc) ha- 
los between hja ^ 0.6 - 0.7 and c/a ~ 0.45 - 0.55. At 
the same time, we find that the halo figure tumbles very 
slowly, ~ 77 over the Hubble time. The DM collapse 
pumps most of the angular m,omentum, into the internal 
circulation within the halo, while the baryonic collapse 
is stopped by the centrifugal barrier in a flat disk. Thus 
the baryons, due to dissipation, collapse much more pro- 
foundly. 

The halo triaxiality decreases with time during the disk 
growth. Two processes contribute to this — (1) the in- 
creasing central mass concentration in the model and the 
appearance of of numerous scattering centers (blobs) in 
the halo, and (2) the out-of phase response of the baryons 
in the disk to the halo shape. The first effect reduces 
somewhat both the fiatness and prolateness of the halo, 
even in the absence of a baryonic component, much more 
so in the presence of baryons. The second effect results 
from the negligible tumbling of the halo figure — under 
these circumstances the inner ILR (if it exists) and the 
outer ILR are pushed to the center and to large radii, 
respectively, and the disk acquires an elongation perpen- 
dicular to the major axis of the halo, thus diluting its 
potential in the equatorial plane. This effect is model 
independent and is unrelated to dissipation. 

The angular momentum (J) flow within the evolving 
models has been always directed outwards — the inner 
disk is losing J to the outer disk and to the surrounding 
DM halo. Overall, the baryons lose their J while the 
DM increases it in the process of evolution. On the av- 
erage, the DM halo acquires about 2.5% -3% of J while 
the baryons lose about 25%-30% of their original mo- 
mentum. The total J is conserved to less than 0.5% 
even in the models with SF. The specific angular mo- 
mentum (j) behaves in a more complicated manner — 
for the DM it is almost constant over the Hubble time. 
For baryons, j decreases with time but each of the com- 
ponents separately, the gas and stars, increase their j 
substantially over the time. This happens because the 
first stars form in the center with a minimal j, and the 
first gas consumed by the SF has small j as well. 

The disk evolution can be roughly divided into two 
stages: in the early stage, the baryons flow towards the 
forming small and amorphous disk along a number of 
radial string patterns — this characteristic configuration 
forms in all models and resembles the Cat's Cradle. It 
survives for about 1 Gyr. In the second stage, the disk 
periodically changes its shape from less to more oval and 
frequently develops a pair of grand-design spiral arms. 
The disk shape and the strength of these arms depends 
on their relative angle with the major axis of the DM 
halo. In nearly all models, the disk grows from inside 
out and is gas-dominated initially. About 1/3 to 1/2 of 
the original baryons remains outside the SF region of the 
disk or in the halo in form of hot gas at the end of the 
simulation (see also Sommer-Larsen 2006). 

The star formation (SF) is typically strongly concen- 



trated toward the disk midplane. The gaseous com- 
ponent in the disk is puffed up by the feedback from 
the stellar evolution. We find that the gas layer nar- 
rows substantially when the SF rate there drops below 
~ 5 M0 yr~^ — which typically happens after the first 
5 Gyr. The SF delineates various morphological features 
like early bars, spiral arms and nuclear rings. Late 
bars appear gas-poor and do not show much SF. 

We find that parameters which characterize the disk 
morphology fit the observations of galactic disks in the 
nearby universe, e.g., their radial and vertical scale- 
lengths, shape parameters, bulge-to-disk ratios, etc. The 
typical radial scalelength in the disk increases by a fac- 
tor of 2 over the Hubble time, ending within the range 
of ~ 2 — 3 kpc. During this time the disk grows by more 
than a factor of 10 in mass. 

A caution should be taken when tackling the evolution 
of a coUisionless component (stars) within the central 
few 100 pc, specifically the formation of bulges in disk 
galaxies. Insufficient timestep resolution can lead to be- 
ing unable to resolve the orbits of individual stars — this 
leads to energy non-conservation and overall heating in 
the region. Such heating has been detected in the central 
few 100 pc of our models. It did not affect the large-scale 
evolution but could, in principle, modify the nuclear bar 
strength as well as increase the dispersion velocities in 
the disky component of the bulge. We, therefore, post- 
pone the discussion of the bulge formation, the bulge-disk 
decomposition and the evolution of characteristic scale- 
lengths in our models to a separate publication where 
high-resolution timestep models are presented. 

Essentially, the disks obtained in our simulations range 
from being bulge-dominated to nearly bulgeless. Mod- 
els with larger feedback from the stellar evolution form 
smaller bulges and vice versa. Models with shorter 
timescale for SF formed 'fatter' disks but not necessary 
more massive bulges. Moreover, reducing the critical 
density for SF has led to a distributed SF activity over 
the inner 15 kpc, as opposed to the inside out disk growth 
in other models, and to 'fat' short (small radial scale- 
length) disks. 

All models are characterized by an extensive bar- 
forming activity within the central few kpc. This is espe- 
cially true for the initial few Gyr of the disk formation. 
The early bar axial ratios vary in tandem with the mutual 
bar-halo orientation. These bars appear to be gas-rich, 
channel their gas contents towards the central kpc and 
weaken substantially over this time. In some cases their 
strength was observed to revive, even after a prolonged 
period of time. In other cases, small nuclear bars de- 
veloped in addition to the primary bars. Dynamically, 
we find that these nuclear bars either corotate or tumble 
differently than the large bars. In all the latter cases, we 
find that this was preceded by a gas inflow to the region 
and that nuclear bars tumble substantially faster than 
their large-scale counterparts. Some of the nuclear bars 
appear not to have visible stellar components — they 
quickly collapse to the center. Nuclear rings and spirals 
appear frequently over the Hubble time and the former 
ones are robust features in the disk. 

Some models show large inflow rates of the gas into 
the central few kpc leading to the Jeans instability in the 
gas and the formation of massive clumps. These clumps 
experience dynamical friction within the disk and spiral 
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in before being consumed by the SF process. Such a 
regime has been investigated by Shlosman & Noguchi 
(1993) who found that the bar growth is suppressed due 
to the disk being heated by the clumps (see also Immeli 
et al. 2004). Here we confirm that indeed the bars have 
been substantially weakened by these clumps, or did not 
form in the extreme cases. 

Our main results indicate that it is possible to form 
disks which are in agreement with those in the local 
universe, e.g., with observed sizes, scalelengths, shapes, 
mass distributions, etc. Furthermore, we obtain rotation 
curves that are flat and with the disk/halo contributions 
as observed (e.g., Bosma, 2004). The so-called angular 
momentum 'catastrophe' discussed in section 1 is natu- 
rally avoided in this approach. The baryons lose only 
about 25% of their angular momentum because both the 
original J present and the feedback from the stellar evo- 
lution reduce their accretion rate onto the disk where 
most of the J is lost. We find that in all models the to- 
tal angular momentum of baryons (stars-|-gas) within the 
disk region closely follows that of the DM there. After 
the initial ^ 3 Gyr, the baryons rapidly lose their J that 
flows to the DM until 'equipartition.' 

Moreover, a full range of bulge-dominated to nearly 
bulgeless disks have been obtained as a result of the evo- 
lution over the Hubble time. The disks formed in our 
numerical simulations do not appear to suffer from be- 
ing too compact and having small radial or large vertical 
scalelengths. The growth of the disk cannot he modeled 
realiahly without resolving its radial scalelength — an is- 
sue which is difficult to address at present in cosmological 
simulations and which is resolved in our approach. The 
problem of over-cooling — when the gas radiates away 
the energy deposited by feedback from the stellar evolu- 
tion (OB stellar radiation-driven winds and SN), is also 
addressed here (see section 2). 

Bars provide the strongest impetus for an internal evo- 
lution in disk galaxies by imposing gravitational torques 
which trigger the mass and angular momentum redistri- 
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butions within the central ~ 10 kpc. While prolate DM 
halos appear to be incompatible with large-scale stellar 
bars (El-Zant & Shlosman 2002; Berentzen et al. 2006), 
the growth of disks acts to wash out the equatorial pro- 
lateness within, at least, the inner halos (Berentzen & 
Shlosman 2006). The present work confirms this trend 
— we observe a decrease in the halo flatness even well 
outside the disk, and elimination of its prolateness within 
the inner (< 20 kpc) halo and reduction in the outer halo. 
The disk shape and morphology are affected by the halo 
shape during the first few Gyr — the disks appear oval 
in contradiction of those observed in the local universe 
(e.g., Rix & Zaritsky 1995). The interaction between the 
ovally-shaped disks and prolate halos provide for the spi- 
ral driving which clearly depends on the mutual disk-halo 
orientations. 

In summary, we have explored different aspects of disk 
galaxy formation in the framework of the collapse of 
an isolated density pertubation embedded in the Hub- 
ble flow. This allows us to focus on specific processes 
normally unresolved in the cosmological simulations. On 
the other hand, we lose the effect of the substructure on 
the model. Specifically, we have aimed at understand- 
ing the imprint of various feedback parameters describ- 
ing the star formation on the basic charateristics of the 
disk-halo system. One should consider this effort as an 
exploratory one, which requires a full implementation of 
realitic initial conditions in tandem with a higher numer- 
ical resolution. 
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